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Abstract. Satellite Synthetic Aperture Radar (SAR) observations can provide a global 
view of ocean swell fields when using a specific "wave mode" sampling. A methodology 
is presented to routinely derive integral properties of the longer wavelength (swell) por- 
tion of the wave spectrum from SAR Level 2 products, and both monitor and predict 
their evolution across ocean basins. SAR-derived estimates of swell height, and energy- 
weighted peak period and direction, are validated against buoy observations, and the peak 
directions are used to project the peak periods in one dimension along the correspond- 
ing great circle route, both forward and back in time, using the peak period group ve- 
locity. The resulting real time dataset of great circle-projected peak periods produces 
two-dimensional maps that can be used to monitor and predict the spatial extent, and 
temporal evolution, of individual ocean swell fields as they propagate from their source 
region to distant coastlines. The methodology is found to be consistent with the disper- 
sive arrival of peak swell periods at a mid-ocean buoy. The simple great circle propa- 
gation method cannot project the swell heights in space like the peak periods, because 
energy evolution along a great circle is a function of the source storm characteristics and 
the unknown swell dissipation rate. A more general geometric optics model is thus pro- 
posed for the far held of the storms. This model is applied here to determine the atten- 
uation over long distances. For one of the largest recorded storms, observations of 15 s 
period swells are consistent with a constant dissipation rate that corresponds to a 3300 km 
e-folding scale for the energy. In this case, swell dissipation is a significant term in the 
wave energy balance at global scales. 



1. Introduction 

Storms over the ocean produce long surface gravity waves 
that propagate as swell out of their generation area. In deep 
water, the wave phase speed C and period T are propor- 
tional. As the phase speed of the dominant waves C p does 
not exceed 1.2 times the wind speed at 10 m height Uio 
[Pierson and Moskowitz, 1964], the longest period waves 
must be generated by very intense winds. For example, the 
generation of waves of period T larger than 16 s requires 
winds with speeds over 18 m s~ x blowing over a distance 
of the order of 1000 km, to produce a significant energy, or 
yet stronger winds over a shorter fetch [Munk et al., 1963]. 
Such a large region of high winds is generally associated with 
a smaller storm center from which the long swells radiate. 
Waves further evolve after their generation with an impor- 
tant transfer of energy towards both high and low frequen- 
cies, due to nonlinear wave-wave interactions. Away from 
that core, nonlinear interactions become negligible [Hassel- 
mann, 1963] and long period swells have been observed to 
propagate over very large distances, up to half-way around 
the globe [Munk et al., 1963], radiating a large amount of 
momentum and energy across ocean basins. This measur- 
able long-distance propagation is made possible by a limited 
loss of energy. 
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The wave field at any time t, latitude <j> and longitude 
A, is described by its spectral densities G as a function of 
frequency / and direction 9. In the limit of geometrical op- 
tics, the spectral density is radiated at the group speed C g 
in the direction of wave propagation, and can be expressed 
as a function of G at any previous time to. Allowing for a 
spatial decay at a rate fi, the spectral energy balance is [e.g. 
Munk et al, 1963] 

G(t, <t>, A, /, 6) = G (t , 0o, Ao, /, Oo) exp Qf -^C B Atj (1) 

In deep water without current, the initial position 0o, Ao 
and direction 9o are given by following the great circle that 
goes through the point of coordinates 0, A with a direction 
—6 over a distance X = (i — to)C g = (t — to)g/(A-Kf). This 
corresponds to a spherical distance a — X/R along the great 
circle, where R is the Earth radius. 

Equation (1) can be used to invert n from wave measure- 
ments. For swell periods shorter than 13 s, Snodgrass et al. 
[1966] have measured an e-folding scale L e — 1/ ^ — 5000 km 
(this number corresponds to a 0.1 dB/degree attenuation in 
their analysis). For larger periods, Snodgrass et al. [1966] 
could only conclude that L e is larger, possibly infinite. In 
the past 40 years, little progress has been made on these con- 
clusions [WISE Group, 2007]. Yet this question if of high 
practical importance, either for wave forecasting [e.g. Ras- 
cle et al. , 2008] or other geophysical investigations regarding 
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air-sea fluxes or microseismic noise [e.g. Grachev and Fairall, 
2001; Kedar et al, 2008]. 

A theoretical upper bound for L e is given by the viscous 
theory [Dore, 1978, see also Appendix A for a simple deriva- 
tion]. According to theory, the largest shears are found right 
above the water surface, and the air viscosity dominates the 
dissipation of swells, giving, in deep water, 



where v a is the air viscosity, a — 2n/T. For T — 13 s this 
gives L e ,max = 45000 km, which means that over a realistic 
propagation distance of 10000 km the energy of 13 s swells 
is only reduced by 25%. 

Swells are thus expected to be very consistent over dis- 
tances that are only limited by the size of ocean basins. The 
analysis of swells at this global scale should provide insights 
into their dynamics, including propagation and dissipation, 
but also into the structure of the generating areas, in a way 
similar to the use of the cosmic microwave background for 
the analysis of the early universe. 

The present paper provides two important intermediate 
steps toward that goal. First, we demonstrate in section 2 
how sparse data from a single space-borne synthetic aperture 
radar can be combined dynamically to provide a consistent 
picture of swell fields. This internal consistency reveals the 
quality of the SAR-derived dataset which we further verify 
quantitatively with buoy data. In section 3, we discuss and 
derive the asymptotic far-field swell energy evolution. Nu- 
merical investigations are performed to check the validity of 
the asymptotic solutions. This result provides a tool to in- 
terpret measured swell heights in terms of propagation and 
dissipation. This method is illustrated with one example 
that corresponds to a strong swell dissipation. Conclusions 
follow in section 4. 

2. Space-time consistency of space-borne 
swell observations 

Investigations by Holt et al. [1998] and Heimbach and 
Hasselmann [2000] have shown that space-borne synthetic 
aperture radar (SAR) data can be used to image the same 
swell field over 3 to 10 days as it propagates along the ocean 
surface. These preliminary studies have shown that the com- 
bination of SAR data at different places and times yields a 
position of the generating storm, and predictions for the 
arrival time of swells with different periods and directions. 
Heimbach and Hasselmann [2000] have further pointed to 
shortcomings in the wind provided by an atmospheric cir- 
culation model for a given Southern Ocean storm, based on 
systematic biases in wind-forced wave model results com- 
pared to SAR observations. Unfortunately, the systematic 
analysis of such data has been very limited, and generally 
confined to data assimilation in wave forecasting models [e.g. 
Hasselmann et al. 1997; Breivik et al. 1998; Aouf et al. 
2006]. This narrow use of SAR data is due to three essential 
difficulties. 

First, a SAR image is not a picture of the ocean surface 
and the relationship between the spectrum of the SAR image 
and that of the ocean surface elevation is nonlinear and fairly 
complex [e.g. Krogstad, 1992]. Sophisticated methods have 
been developed in order to estimate the surface elevation 
spectrum [e.g. Hasselmann et al., 1996; Schulz-Stellenfleth 
et al., 2005]. These methods had to be implemented by 
the user of the data, and generally required some a priori 
first guess of the wave field provided by a numerical model. 
For longer wave systems, the imaging mechanisms are es- 
sentially quasi-linear, making possible a simpler methodol- 
ogy used by the European Space Agency (ESA) to generate 
a level 2 (L2) product. The method is fully described by 



Chapron et al. [2001]. It uses no outside wave information, 
and builds on the use of complex SAR data developed by En- 
gen and Johnsen [1995] to remove the 180° directional am- 
biguity in wave propagation direction. The quality of the L2 
data has been repeatedly analyzed [e.g. Johnsen et al, 2006; 
Collard et al, 2005]. Because long ocean swells have large 
wavelengths and smaller steepnesses, the L2 products corre- 
sponding to this spectral range have higher relative quality, 
confirming that the imaging mechanism is well described 
under the quasi-linear assumption. 

All SAR data used here are such L2 products, provided 
by ESA and obtained with the L2 processor version op- 
erational since November 2007, and described by Johnsen 
and Collard [2004]. The data for times before that date 
were reprocessed with this same processor. In previous real- 
time data, frequent low wavenumber artefacts were caused 
by insufficient filtering of non-wave signatures in the radar 
images. This filtering is necessary to remove the contribu- 
tions of atmospheric patterns or other surface phenomena 
like ships, slicks, sea ice, or islands, with spectral signa- 
tures that can overlap the swell spectra. The L2 product 
contains directional wave spectra with a resolution of 10° in 
directions and an exponential discretization in wavenumbers 
spanning wavelengths of 30 to 800 m with 24 exponentially 
spaced wavenumbers, corresponding to wave periods with a 
7% increment from one to the next. 

The second practical problem is that the data obtained 
from an orbiting platform are sparse and with a sampling 
that makes a direct analysis difficult. Hereafter, we show 
that the space-time consistency of the swell field can be used 
to fill in the gaps in the observations and produce continu- 
ous observations of swell periods and directions in space and 
time. 

Third, and last, for a simple use of SAR data, some pa- 
rameters that are not affected by the variable resolution is 
SAR scenes [Kerbaol et al., 1998]. Schulz-Stellenfleth et al. 
[2007] have proposed to produce parameters representing the 
entire sea state by extending the resolved spectrum with an 
empirical windsa contribution. Here we take the opposite 
approach and restrict the resolved part of the spectrum by 
using a spectral partitioning (see Appendix B for details) to 
retrieve the swell significant wave height H as , defined as four 
times the square root of the energy of one swell system, and 
the peak period T p and peak directions 8 P . Thus one SAR 
typically produced one or two swell parameters for distinct 
swell systems. 

For H as , only very limited validations have been per- 
formed [Collard et al, 2006]. We thus perform a thorough 
analysis of SAR-derived H ss comparing to co-located buoy 
data (see Appendix B for details). 

The bias on H ss derived from ESA Level 2 products is 
found to be primarily a function of the swell height and 
wind speed, increasing with height and decreasing with wind 
speed. Variations in standard deviation are dominated by 
the swell height and peak period, with the most accurate 
estimations for intermediate periods of 14 to 17 s.For wind 
speeds in the range 3 to 8 m s _1 , H ss has a bias of 0.24 m 
and the standard deviation of the errors is 0.29 m. 

We thus corrected the values of H ss by subtracting a bias 
model given by 

b h = 0.11 + 0.1H SS - 0.1 max{0, Uiosar - 7} (3) 

where H as is in meters and the wind speed Uiosar is in 
m s _1 . From now on, all the reported values of H ss will be 
corrected using this expression. After correction, the stan- 
dard deviation of H as estimates is reduced to less than 

o h = 0.10 + min{0.25# ss , 0.8} (4) 
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where ah and H sa are in meters. 

The quality of T p and 9 P has already been carefully stud- 
ied and are routinely monitored. The root mean square 
(r.m.s.) error on T p is less than 10% of the measured value 
for T p > 12 s, and the r.m.s. error on 9 P is 22° for these 
same swells [Johnsen and Collard, 2004], with little bias. 
Because few directional in situ measurements are available, 
we demonstrate here an original semi-quantitative dynami- 
cal validation of these two parameters. 

2.1. Virtual wave observers 

Given these SAR-derived estimates of 9 P and T p , linear 
dispersion relationship and the principle of geometrical op- 
tics can then be exploited to predict arrival times and loca- 
tions of the swell. 




interrogation window 

1(2,1) 

virtual buoy 

position 




f 



t (4) 



t 

t (3) 



t (1)/tld.D 



ti(2,2) 



Figure 1. Schematic definition of a virtual wave ob- 
server. Any SAR observation i is available at a time 
to(i) on the black dots. All swell partitions (i, j) (here 
indicated by the arrows) are propagated and may cross 
the interrogation window from time t\(i,j) to £2(2, j). 



In order to obtain swell conditions at the location of a 
"virtual wave observer", we define an interrogation window 
covering 2 by 2 degrees in latitude and longitude. Accord- 
ing to the SAR-derived peak parameters, swell partitions 
from the entire ocean basin are propagated, both forward 
and backward, along great circles in space and time. These 
theoretical trajectories are followed with a constant group 
speed gT p /(4:ir), starting off in direction 9 P from the obser- 
vation point. From any observation time to, these great cir- 
cles may cut through the interrogation window from times 
ti to t2 (thick solid lines in figure 1, with different colors 
for different partitions). As the maximum value of |t2 — to 
and \ti — to I is increased from 6 hours to 6 days, the time- 
evolution of the peak frequencies and peak directions at the 
virtual observer gradually reveals similar ridges to the one 
observed in real buoy measurements (figure 2). 

In fig. 2b-d, each horizontal colored segment corresponds 
to one swell partition that crosses the spatial window be- 
tween times ti and t^- Some segments are very short, cor- 
responding to trajectories that barely cut one corner of the 
window. 

Clearly the SAR detects the direction of the most en- 
ergetic part of the wave spectrum measured by the buoys 
(fig. 2a). At frequencies above 0.1Hz, the virtual observer 
patterns appear rather noisy. Shorter scales are not so cor- 
related. These shorter components are often observed as 
part of the wind sea for which the propagation with a sin- 
gle group speed and direction is not a good approximation. 
Also, propagated high frequency swells, such as the 0.12 Hz 



waves coming from direction 200 on July 10, do not show 
up in the real buoy record. This is possibly the result of a 
relatively high dissipation rate for these swells. 
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Figure 2. (a) Energy and mean direction spectrum mea- 
sured in situ by the Christmas Island buoy (WMO num- 
ber 51028): contours, equally spaced from 0.1 to 1.4, indi- 
cate the natural logarithm of the spectral energy density 
F(f). Colors indicate the mean arrival direction at each 
frequency, (b) to (e) Peak direction (colors) as a function 
of time and peak frequency for swell partitions at a SAR 
virtual buoy located around the Christmas Island buoy 
(WMO number 51028). The maximum propagation time 
to produce the virtual buoy data is increased from 6 hours 
(b) to 6 days (e). The sloping straight line fitted to the 
observed SAR ridge from July 16 at 0.05 Hz to July 21 at 
0.105 Hz is the same as line in (a) that corresponds to the 
buoy observation. This delayed arrival would correspond 
to a point source at 6100 km from Christmas Island. 



For frequencies below 0.08 Hz, the virtual observer shows 
ridge-like structures similar to those observed in situ due 
to the dispersive arrivals of swells from remote storms [e.g. 
Munk et at, 1963; Gjevik et at, 1988]. Even the faintest 
events are well detected, such as the 0.06 Hz arrival on 23 
July, even though that 15 s swell of 0.5 m is dwarfed by a 
another 0.8 m and 12 s swell, and a 2 m and 8 s wind sea. 
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Swell detection with the virtual observer reaches its limits 
when the swell height is very low, such as on 3 July with a 
20 s 0.3 m swell well detected by the buoy (the green-orange 
ridge at 0.05 Hz). Because we use a single trajectory ema- 
nating from one observed swell partition the relatively small 
interrogation window can easily be missed after 10000 km 
of propagation. 

Our discrete propagation technique suffers from a ran- 
domized version of the "garden sprinkler effect" that, if not 
corrected for, can create unrealistic flower-like patterns in 
the far field of storms in numerical wave models that use 
a discretized spectrum [e.g. Tolman, 2002]. Our choice of 
a single group speed and direction, because a narrow swell 
spectrum is not resolved by the SAR, produces a discrete 
wave field (the dots in figure 3). With the present process- 
ing this is smoothed by the finite size of our interrogation 
window (figure 1). An extension of the present technique 
could use neighboring group speeds and directions to take 
into account the frequency and directional spread of the sea 
state, which would allow the use of a smaller window. Just 
like the estimation of propagated wave heights, discussed 
below, the estimation of a spectra width that cannot be re- 
solved by the SAR may use some further information on the 
structure of the generating storm. 

Other errors in the present technique can also be at- 
tributed to the SAR processing. In particular, a maximum 
value is defined for the transfer function used to obtain the 
wave spectrum from the SAR image [Johnsen and Collard, 
2004]. Although this is designed to prevent the amplifica- 
tion of measurement noise, long swells such as this 20 s event 
have very small slopes, and it is likely that they are under- 
estimated in the wave spectrum due to this threshold in the 
processing. 

When propagated for 6 days, without any other informa- 
tion than the peak frequency and direction at the time of 
observations, the waves are remarkably consistent with the 
latest local observations. For the southern swells arriving 
at Christmas Island between July 16 and 21 (figure 2), the 
difference in arrival times given by the virtual oberver and 
real buoy is typically less than 12 h. This is less than 10% 
of the maximum time between the SAR observation and the 
virtual observer record. This implies that the accuracy of 
the peak period estimate for each SAR partition must also 
be less than 10%, consistent with previous validation stud- 
ies [Johnsen and Collard, 2004; Johnsen et al., 2006]. The 
consistency of the arrival directions along the ridges also 
suggests that the root mean square (RMS) error in peak di- 
rection estimates must be close to 20°, comparable to the 
22° RMS difference between mean wave directions obtained 
from SAR wave mode and a numerical wave model for waves 
with periods longer than 12 s [Johnsen and Collard, 2004]. 

Although it cannot replace the spectral resolving power of 
a buoy, the performance of the virtual observer is therefore 
comparable or better to that of human observers in terms of 
peak period and direction [Munk and Traylor, 1947]. The 
really missing bit is a wave height estimate along the swell 
propagation path. We will show that this may be obtained 
by estimating the source storm characteristics and the dis- 
sipation rate of swell energy. 

2.2. Storm source identification and "fireworks" 

Along the estimated trajectories, virtual observations can 
further be produced in a similar fashion. The animation of 



these propagated swells confirms the very well organized na- 
ture of storm swells crossing large ocean basins. 
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Figure 3. Snapshot of the "fireworks" animation given 
in the auxiliary material, for July 17 2004 at 0:00 UTC. 
Each of the 1071 colored dots represent one observed 
swell partition, within 6 days of its observation, displaced 
along a great circle with the group speed corresponding 
to the detected peak period in the direction of the de- 
tected peak direction. Only swells with tracks that passes 
within 1000 km of the storm center (red disk) have been 
retained. 

From the relatively sparse and track-based initial satel- 
lite observation sampling, the swell persistency can then be 
used to capture "fireworks" patterns exploding from the few 
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Figure 4. Finding the source storm. All swells with a 
17 ± 0.5 s period that were identified in 13 days of EN- 
VISAT synthetic aperture radar data over the Pacific, 
are re-focussed from their location of observation (filled 
dots) following their direction of arrival at the theoreti- 
cal group speed for 17 s waves. This focussing reveals a 
single swell generation event, well defined in space and 
time (pink to red disks). The back-tracking trajectories 
are color-dated from black (July 9 2004 18:00 UTC) to 
red (July 22 2004 18:00 UTC). ' 
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intense storms that occur over a period of several days (see 
figure 3 and auxiliary animation). In large ocean basins 
where swells are likely to be imaged several times by the 
same satellite, these fireworks can be used to estimate the 
time of arrival of swells from any given storm [e.g. Holt 
et al, 1998]. For this reason, these animations have been 
produced routinely every day since August 2007, for the 
Pacific, Indian, and Atlantic oceans (see http://www. boost- 
technologies, com/esa/images/, e.g. nrt_pac.gif for the Pa- 
cific Ocean). 

Using backward trajectories, the location and date of 
swell sources can further be defined as the spatial and tem- 
poral center of the convergence area and time of the trajec- 
tories. These positions have been verified to correspond to 
high wind conditions observed by scatterometers and repro- 
duced by ECMWF wind analyses. We consider these storms 
to be the source of all the swell partitions that produce tra- 
jectories that pass within 12 hours and 2000 km of their 
center. This processing, similar to the one performed by 
Heimbach and Hasselmann [2000] , provides a global view of 
swell fields in both space and time, extending the coverage of 
similar techniques based on buoy data [Hanson and Phillips, 
2001]. In figure 4, a swell covers one Earth quadrant away 
from the storm, with a large detection gap that extends from 
the Southern Pacific to California. This blank area is the 
long shadow cast by French Polynesia where wave energy is 
dissipated in the surf [e.g. Munk et al., 1963]. Observations 
were restricted to swell partitions with periods close to 17 s, 
but the full dataset typically covers swells with periods of 
12 to 18 s, as shown in figure 3. 

The apparent self-consistency of both the virtual buoy 
plot (figure 2.d) and the fireworks animations, are the result 
of the large auto-correlation length of the swell fields, which 
was expected from the in situ measurements of Darbyshire 
[1958]; Munk et al. [1963]; Snodgrass et al. [1966]. Yet, these 
plots could not exhibit such patterns without a good accu- 
racy of the SAR-derived peak periods and directions, used 
in the propagation methodology. 

3. Far-field swell energy 

All these illustrations of forward-backward ray tracing in- 
dicate the potential to use a simple Geometrical Optics (GO) 
strategy. The next goal is then to determine the strength 
of the far-field radiated swell energy. This requires the def- 
inition of a swell source, and an estimation of the swell en- 
ergy dissipation scale. For this we define the time to as 
an initial condition after which there is no significant wave 
generation or non-linear evolution, for frequencies less than 
/max- Namely, at to all the wave components with smaller 
frequencies have already been generated, so that the radi- 
ation of these waves is essentially fossil and fully governed 
by geometrical optics. The possible effects of diffraction and 
scattering are discussed by Munk et al. [1963], and, together 
with dissipation, will cause deviations from the G.O. model 
outlined below. We therefore make no restrictive hypothe- 
sis on what happens before to, and thus the motion of the 
generating storm has no direct effect on our results, but it 
obviously modifies the spatial distribution of the energy at 
to, which will be relevant. 

In reality, swells evolve over the course of their propaga- 
tion as the result of their interactions with the local winds, 
mutual wave- wave interactions, interactions with other wave 
systems, including the local wind sea. Swells are also ex- 
pected to evolve according to interactions with other oceanic 
motions that affect the upper ocean, namely surface cur- 
rents, internal waves [e.g. Kudryavtsev, 1994] and turbu- 
lence. Depth and island scattering effects must also care- 
fully be taken into account [Snodgrass et al, 1966; WISE 



Group, 2007]. Compared to these different mechanisms, fre- 
quency dispersion and angular spreading effects are certainly 
the first leading order phenomena to take into account for 
the major part of the decrease in the height of the swell 
systems. Indeed, as the long swell systems will be charac- 
terized by relatively small steepness parameters, nonlinear 
mutual wave-wave interactions do not appear to be impor- 
tant in scattering surface wave energy more than a few storm 
radius distances outside an active generating area Hassel- 
mann [1963]. Furthermore, the level of turbulence in the 
ocean does not appear to significantly affect the waves Fab- 
rikant and Raevsky [1994]; Ardhuin and Jenkins [2006], and 
the conversion of surface wave energy into internal gravity 
wave energy by wave-wave interactions does not seem to 
be a leading order sink term for the energy balance of sur- 
face gravity waves. Finally, for the very long trans-ocean 
fast propagating swell components, surface current bending 
effects, proportional to the ratio between vertical current 
vorticity and the group velocity, may also be considered as 
residual effects. Away from island obstructions, the ratio 
between the angular width along the great-circle observa- 
tory points at very large distances from the generating area, 
and the mean spread in the generating area is approximately 
proportional to 1/ sin a with a the spherical distance from 
the storm, center The change in spectral density F defined 
by 

/•2tt 

F(t,4>,\,f) = / G(t,4>,\,f,0)d0 (5) 
Jo 

follows the spatial expansion (close to the source) and con- 
traction (as waves approach the antipodes for a > tt/2) of 
the energy front. This transversal dispersion is associated 
to a narrowing of the directional spectrum G, for a < 7r/2, 
and a broadening for larger distances. 

This approximation applies to large distances, and rela- 
tively small source regions. Closer to the source, the ap- 
proximation does not hold. Swell amplitudes radiating from 
large extended sources will decrease more slowly than swell 
amplitudes emanating from compact sources. 

Moreover, we can represent swell waves as a linear su- 
perposition of harmonic waves in narrow spectral band. 
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Figure 5. Dispersion of linear waves in one dimension. 
At any given time the spectrum is given by a propagation 
of the spectra at t = to. Taking xi = xo + (ti — to)C g (f p ) 
the spectral density at f p is the same as for time to, but 
the spectrum is narrower which gives a smaller elevation 
variance, E(ti,xi) < E(to,xo). 
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Quite naturally, through the method of stationary phase, the 
group velocity is defined and the slowly-varying wave enve- 
lope is found to decay. This decay is inversely proportional 
to the square-root of the distance (figure 5). Accordingly, 
far away from the generating sources, and in the absence of 
dissipation, the swell energy 



E 3 {t,cj>,\) 



F(t,4>,\,f)df. 



(6) 



should decrease asymptotically like l/(asina) when follow- 
ing a wave group (see Appendix C for a detailed proof). 

3.1. The Snodgrass et al. [1966] method 

Using measurements with a limited or no directional res- 
olution, Snodgrass et al. [1966] assumed that wave propa- 
gation was completely blocked by waters shallower than 60 
fathoms (approximately 110 m), and that diffraction could 
be neglected. For example, in figure 6, the island would 
be represented by the 60 fathom depth contour. These au- 
thors then estimated a loss of swell energy from the de- 
viation of the ratio of directionally-integrated spectra, e.g. 
F(x2, fo)/F(xi, fa), compared to what is expected from the 
lateral dispersion effect, taking into account islands. Cur- 
rents, shallow water areas and diffraction effects around is- 
lands are neglected here. These effects are discussed by 
Snodgrass et al. [1966]. 
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Figure 6. Dispersion of linear waves in two dimension, 
represented here on a flat surface for simplicity. The vari- 
able g has been dropped as the spectra shown here are 
all at g — 0. We call "source" the region of the ocean 
where waves with frequencies smaller than / m i n can be 
found. As time goes by, the source expands in space due 
to both frequency dispersion (like in ID), and geometrical 
dispersion. The wave energy with frequency /o that will 
be observed at point xi (respectively xz) at time t\ (re- 
spectively £2) is, at time to, along the thin arc circle AB 
(respectively the thick dotted arc circle CD). Due to the 
small island between X\ and X2, the energy that would 
have been recorded at x%, if the island were not present, is 
actually dissipated on the shore of the island. As a result 
the local energy density E(x2) is reduced. At frequency 
fo, contributions to E{x2) only come from angles #5 to 
f?e. The directional spectra (bottom) are thus affected 
by the blocking effect of islands, and the directional nar- 
rowing as one goes further from xo (on the Earth this 
narrowing reverses after 10000 km of propagation, due to 
the sphericity). 



Rigorously, their method is inexact because the record- 
ing stations Xi and x% measure wave groups that had neither 
exactly same propagation directions nor the same position 
when they were near xq. Yet, because the wave field in 
the neighborhood of xq is the superposition of many in- 
dependent wave trains, one can assume that the spectral 
density F is a smooth function of the direction. Then we 
may say that over the intervals #3 to O4 or 85 to 6*6 , G does 
not change so much, i.e. in figure 6, G(to, cj>n, Ad, /, 9e) ~ 
G(to, 4>d, Ad, /, 9$). On the sphere the application of eqs. 
(1) and (5) yields 



F{t,4>,\fo) 



1 



G(t ,4>,\, fo,0o) 



l + O 



\ sm a J 



(7) 



where the integral is performed over the line segment joining 
C to D. The error relative to the asymptote is the sum of 
two terms. One is proportional to the spatial gradient of F, 
due to the change from the arc circle to the segment, and the 
other corresponds to the relative variation of G with 9 over 
the range #5 to 9e, which is small in the far field, provided 
that the directional spectrum is smooth enough. 

Under these two smoothness assumptions, and for large 
sin(a), the ratios of spectral densities F at xi and x%, as 
used by Snodgrass et al. [1966], can be compared to the con- 
servative equations (7) and used to diagnose the dissipation 
of swell energy. This does not apply if x\ is in the vicin- 
ity of the storm or its antipode, where the observed arrival 
direction span a large range. 

In practice, this method can be very sensitive to the cor- 
rect estimation of the island shadowing. For that reason, the 




N 



Figure 7. Geometry of the "fossil' swell field distribu- 
tion at time to = (shaded area) and observation condi- 
tions. Any point P of colatitude <p and longitude A inside 
of the storm, generates waves that are observed at point 
O. At time t the observed waves that come from P have 
a well defined frequency given by eq. (C2), function of 
the spherical distance a' between P and O, and a well de- 
fined direction 9 at P, relative to the North, which gives 
a direction tt — 9' at O. In the triangle OPS the angles 
A, 9' and n — 9 are related to the distances a' , f> and a by 
the usual spherical trigonometry relationships, e.g. eq. 
(C4). 
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measurement route chosen by Snodgrass et al. [1966] was far 
from ideal. Because they needed land to install most record- 
ing stations for the wave measurements, they used an almost 
north-south great circle that extends from the south of New- 
Zealand (Cape Palisader) to Alaska (Yakutat), a route pep- 
pered with islands in its southern part, and partially blocked 
by the Hawaiian chain in its northern part. Also, storms typ- 
ically refuse to line up with any measurement array. Their 
pre-defined great circle, although designed to follow a typi- 
cal Southern winter swell propagation path, always deviated 
by some extent form the actual track followed by the most 
energetic swells they recorded. For the Indian ocean storms, 
this difficulty was reduced by the relatively narrow range of 
angles that allows propagation from the Indian to the Pacific 
ocean. 

3.2. A method using global swell heights 

Now using an instrument with a global coverage, we 
can carefully avoid both problems by choosing propagation 
paths far away from the smallest island, and by exploiting 
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Figure 8. (a) Location of SAR observations with a 15 s 
peak period swell system corresponding to the 12 Febru- 
ary source, with outgoing directions of 74 to 90°. The 
same swell was also observed at all buoys from 46075 off 
Western Alaska, to 51001 in Hawaii. The dash-dotted 
line represent great circles leaving the storm source with 
directions 42, 59, 74, 90 and 106°. (b) Observed swell 
wave height as a function of distance. The solid lines 
represent theoretical decays using no dissipation (blue), 
or a the fitted linear dissipation (green), for swells ob- 
served in February 2007. Circled dots are the observa- 
tions used in the fitting procedure. Error bars show one 
standard deviation of the expected error on each SAR 
measurement. 



only observations well aligned with the storms. However, 
due to the limited spectral resolution inherent to the SAR 
wave mode image size and processing, we cannot use the 
spectral distribution G or F of the energy, and can only use 
the energy E s integrated over a swell partition. We thus 
need for E s the equivalent of eq. (7) for F. 

For simplicity, and without loss of generality, we take the 
source storm centered at time to = on the pole S defined 
by a colatitude ip = (f> — n/2 — 0, so that the distance from 
the storm center is r = Rip. We consider the swell energy E s 
observed at a position defined by the spherical distance a 
and we take the reference meridian to be in the direction of 
the observation point (figure 7). We will later assume that 
the source area is relatively small with a size RA a , where 
A a is the maximum value taken by ip (figure 7). In all the 
following derivations, we have chosen a fixed frequency fo 
and we follow a wave group of that frequency. The time of 
observation t is thus related to a by t — Ra/Cg(fo), so that 
the variable t will be omitted. 

In appendix C, we prove and verify that, in the absence 
of dissipation the swell energy E 3 decreases like l/(asin«) 
for large values of a, and that typical storms should pro- 
duce swells within 20% of this asymptote at distances larger 
than 4000 km from the storm center. A much larger ob- 
served deviation should thus reflect a gain or loss of energy 
by the propagating swells. The expected departures from 
the asymptotic evolution should be compared to those due 
to swell dissipation or generation. Even with perfect SAR 
observations, this is the intrinsic limit of the present method. 
A 20% error in energy conserving conditions may be misin- 
terpreted as a dissipation or generation with an e-folding 
scale of the order of 20000 km, which gives a 20% energy 
change as waves propagate from 4000 to 8000 km away from 
the storm source. 

3.3. Illustration 

To illustrate the method described above, we analyse of 
one of the most powerful swell field recorded over the past 
4 years by ENVISAT's ASAR. The swell case illustrated in 
figure 4 is not well suited due to the islands in the south- 
north swell tracks and the poor sampling of ENVISAT for 
the tracks going north-east from the storm, we have thus 
chosen another source, found on February 12 2007 at 18:00 
UTC, and located at 168 E and 38 N. This swell was gen- 
erated by a storm moving Eastward with maximum West- 
erly winds of 26 m s _1 at the indicated date, and subsiding 
to less that 22 m s -1 and veering to South- Westerly 12 h 
later (according to ECMWF 0.5° resolution analyses). The 
storm motion in the direction of wave propagation certainly 
helped to amplify the local windsea, with a maximum sig- 
nificant wave height H s of 14.1 m at 00:00 UTC on February 
13 (according to a numerical wave model configuration that 
is otherwise verified to produce root mean square errors less 
than 9% for H s > 8 m). 

Using SAR-measured wave periods and directions at dif- 
ferent times and locations, we follow great circle trajectories 
backwards at the theoretical group velocity. The location 
and date of the swell source is defined as the spatial and 
temporal center of the convergence area and time of the 
trajectories. We chose a central peak period, here 15 s, and 
track the swells forward in space and time, starting from the 
source center at an angle 9q, following ideal geodesic paths 
in search of SAR observations. Along each track, SAR data 
are selected if they are acquired within 3 hours and 100 km 
from the theoretical time and position. Great circle tracks 
are traced from the source in all directions, except for an- 
gular sectors with islands. 

In order to obtain enough SAR data, we repeat this op- 
eration for regularly spaced values of 8o with a step of 2°. 
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In our case, when varying 9o from 74 to 90° (counted clock- 
wise from North), this procedure produced 58 SAR measure- 
ments with one swell partition that had a peak wavelength 
and direction within 50 m and 20° of the expected theoret- 
ical value. 

If no energy is lost by the wave field, E s decreases asymp- 
totically as l/[asin(a)] away from the source. Among the 
58 swell observations, we further removed all the data within 
4000 km of the source center, to make sure that the remain- 
ing data are in the far field of the storm, and data with a 
significant swell height H ss less than 0.5 m, after bias cor- 
rection based on the error model. This makes sure that the 
signal to noise ratio in the image is large enough so that the 
wave height estimation is accurate enough. 

We then have 35 observations for which we assume that 
E s is only a function of a, and we define the dissipation rate 

1 d(E s asina) 
» = -E-R Ta • (8) 

Positive values of /i correspond to losses of wave energy (fig- 
ure 8). We then fit an analytical function H ss (a) to the 
data, defined by a constant fj, and H ss (a = 7r/5), i.e. the 
swell height at 4000 km from the source. Here the couple 
H S3 (a = 7r/5) = 4.4 m and [i = 3.7 x 10~ 7 m _1 gives the 
least square difference between the decay with constant /x 
and the observed swell decays. Further, the uncertainty of 
that dissipation rate may be estimated from the known un- 
certainty of the SAR measurement of H 3 , given by eq. (3)- 
(4). A more simple error model, with larger errors based 
on the H a \2 analysis by Johnsen et al. [2006], does not sig- 
nificantly alter this analysis. Using that error model and 
neglecting other sources of error in the present analysis, we 
perturbed the observed swell heights independently to pro- 
duce a 400 ensemble of synthetic data sets. Taking the 16% 
and 84% levels in the estimation of fi, that would correspond 
to one standard deviation for a Gaussian distribution, we 
find that 3.1 x 10~ 7 < pt < 4.0 x 10~ 7 at the 68% confidence 
level. This is the first ever estimation of the uncertainty on 
an observed swell dissipation rate. These values of fi are 
more than twice larger than reported by Snodgrass et al. 
[1966] for smaller amplitude swells. 

The formidable height of 4.4 m at a distance of 4000 km 
was observed by the SAR for all outgoing directions from at 
least 74 to 106° . This same swell was also recorded by buoys 
in the North-East of Hawaii (NDBC buoy 51001), also with 
a peak period of 15 s, and a height of 3.4 m on 16 February 
2007 at 0:00 UTC. That buoy is located 3300 km away from 
the center and in a direction close to 112°. Looking in the 
North-East quadrant of the storm, one also finds a trace of 
the swell at buoy 46005, off the Washington coast (4900 km 
in direction 59°). There the swell was also observed with 
a 15 s period and a maximum height of 3.2 m on February 
17 at 17:00 UTC, similar to the SAR observations for the 
same distance. For directions closer to Northbound, either 
the generation was weaker or the Alaskan islands sheltered 
the coastal buoys. For example, the same swell was also 
recorded by NDBC buoy number 46075, off Shumagin Is- 
land, Alaska, at a distance of 3000 km from the source, in 
the direction 42°. At that buoy, the peak period was 15.0 s 
with a maximum swell height of 1.3 m on 15 February 2008 
at 18:30 UTC. 

Thus the power radiated by the storm is of the order of is 
0.5 TW at 4000 km from the storm center, spread over a 50° 
angular sector. This power is about 16% of the estimated 
3.2 TW annual mean flux that reaches the world's coastlines 
[Rascle et al, 2008]. However, the observed dissipation rate 
corresponds to an e-folding scale of 3300 km for the energy. 
Taking an average propagation distance of 8000 km, only 
160 GW would make it to the shore. If the same dissipation 
rate prevailed closer to the source, then the power radiated 
at 1000 km form the storm center was 1.4 TW. 



We thus expect that the far field dissipation of swells, in 
spite of the small steepness of these swells, plays a significant 
role in the air-sea energy balance. This effect probably ex- 
plains the systematic positive bias for predicted wave heights 
in wave models that neither account for swell dissipation nor 
assimilate wave measurements [see e.g. Rascle et al, 2008]. 

4. Conclusions 

Taking advantages of the satellite observations of un- 
precedented coverage and quality, investigations can repeat 
and complement the pioneering analysis of swell evolution 
performed in the 1960s. Severe storms can generate rela- 
tively broad spectra of large surface waves. But rapidly, 
the redistribution of energy, through linear dispersion and 
nonlinear interaction mechanisms, becomes very effective. 
The initial wind waves become swells outrunning the wind, 
leading to the apparition long-crested systems. The propa- 
gation properties of these surface gravity waves have been 
found to closely follow principles of geometrical optics. The 
consistent patterns of swell fronts dispersing over thousands 
of kilometers was shown to be useful to provide time series 
at "virtual wave observing stations", filling gaps in space 
and time in between the orbit cycles of observation. When 
compared to buoy measurements, the present results give 
an explicit dynamical validation of the SAR-derived spec- 
tral parameters. As the speed of waves in deep water is 
proportional to their period or wavelength, information car- 
ried by the SAR-derived period and direction distributions, 
observed at a fairly large distance from the generating area, 
pertains to the wind conditions existent up to 15 days be- 
fore. 

We also discussed how the swell energy should, in the 
absence of dissipation, decay in the far field of the storm 
like l/(asina) where a is the spherical distance between 
the storm center and observation point. Exploiting that 
property allowed us to estimate a dissipation rate fi of swell 
energy with unprecedented accuracy, establishing that swell 
dissipation can be a significant term in the global wave en- 
ergy budget. 

The proposed methodolgy performed here requires data 
far enough from the source, typically more than 4000 km, 
in order to approach this simple asymptotic behavior. At 
the same time, the swell amplitude should be large enough 
to be accurately measured by the SAR. Some knowledge 
of the spectral shape and its spatial distribution inside the 
storm may be useful to provide better estimates of fj, for 
low dissipation cases, or closer to the storm centers. These 
further analyses will likely benefit from the joint use of data 
from altimeters, SARs, and other sources of spectral wave 
information. 

A more systematic analysis and interpretation of this dis- 
sipation will be reported elsewhere [Ardhuin et al, 2009a], 
with applications to wave forecasting models [Ardhuin et al., 
2008, 2009b]. The parameterization of the dissipation rate 
could also be used to produce a data-based forecasting sys- 
tem, extending our virtual buoy technique to the estimation 
of swell heights, with a forward propagation of observations. 

Going in the opposite direction, toward the storm source, 
it is possible that the analysis of swell fields could provide 
a "new" way of looking into the poorly observed structure 
of severe storms. Because the usual remote sensing tech- 
niques for estimating wind fields either do not work for very 
high winds or are not well validated [e.g. Quilfen et al., 
2006, 2007], the use of far-field swell information may pro- 
vide an interesting complement to the local wind speeds 
and wave heights. The is inverse problem has already been 
formulated by Munk et al. [1963] who already proposed an 
elegant heterodyning technique to push the spatial resolu- 
tion for the estimation of storm location from swell data, 
while Heimbach and Hasselmann [2000] have proposed to 
use wave models to correct wind field errors. The quality 
of the SAR-derived swell parameters that are coming out of 
today's ENVISAT and tomorrow's Sentinel-1, together with 
a good understanding of the swell energy budget, including 
its dissipation revealed here, may finally enable this vision. 
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Appendix A: 
interaction 



Viscous theory for air-sea 



This gives the low frequency asymptote to the viscous decay 
coefficient [Dore, 1978], 



For the sake of simplicity we will consider here the case of 
monochromatic waves propagating in the x direction only, 
and we will neglect the curvature of the surface. For the 
small steepness swells considered here that latter approxima- 
tion is well founded and a more complete analysis is given 
by Kudryavtsev et Makin (2004). For deep water waves, 
the free stream velocity above the waves, just outside of the 
boundary layer is u+(x,t) = — aacos(kx — at), where a is 
the swell amplitude and a — 2n/T is the radian frequency. 
The sub-surface velocity is U-(x,t) = aacos(kx — at) (figure 
9). Due to the oscillations that propagate at the phase ve- 
locity C, the horizontal advection of any quantity X by the 
flow velocity u, given by udX/dx, can be neglected com- 
pared to its rate of change in time dX/dt since the for- 
mer is a factor u/C smaller than the latter, which is typi- 
cally less that 0.1 for the swells considered here. Defining 
u(x, z, t) — {u(x, z, t)) — U— (x, t), where the brackets denote 
an average over flow realizations for a given wave phase. The 
horizontal momentum equation is thus approximated by, 



du 

at 



j_dp 

Pa dx 



du- 



+ G 



(Al) 



where G represents the divergence of the vertical viscous 
and turbulent fluxes of horizontal momentum, 



d 2 u diu'w') 
az A az 



(A2) 



Because the boundary layer thickness 8 is small compared 
to the wavelength, the pressure gradient in the boundary 
layer is given by the pressure gradient above the boundary 
layer, in balance with the horizontal acceleration. This is 
another way to write Bernoulli's equation (e.g. Mei 1989), 

—dp/dx/p a = — a 2 as'm(kx — at) = du+/dt. (A3) 



This yields 



du _ du + 



+ G 



(A4) 



with the boundary condition for z 5, u goes to 2u+(x,t). 
The equation for the horizontal momentum is thus exactly 
identical to the one for the oscillatory boundary layer over a 
fixed bottom with wave of the same period but with an am- 
plitude twice as large. In the viscous case, G = vd 2 u/dz 2 
and one recovers, after some straightforward algebra, the 
known viscous result, i.e., for z > £, 

, z,t) = 2aa [e~ z+ cos (kx — at — z+) — cos (kx — at)] + 0(p a 

(A5) 



a 2 p 
gC g p w 



(A8) 



This result was previously obtained using a Lagrangian ap- 
proach without all the above simplifying assumptions [We- 
ber and F0rland, 1990, their equation 7.1 in which our p 
correspond to their 2a]. The full viscous result is obtained 
by also considering the water viscosity v w , which gives the 
0(pa/p-w) correction for the motion in the air, and the classi- 
cal dissipation term with a decay p vw = —4:k 2 u w /C g , which 
dominates for the short gravity waves. 

The total disspation rate is simply p = p v + p vw . For 
a clean surface with v w = 3 x 10~ 6 m 2 s" 1 , v — 1.4 x 
10 -5 m 2 s _1 , and p a /p w = 0.0013, the two terms are equal 
for waves with a period T = 1.3 s and a wavelength of 2.6 m. 
The air viscosity dominates for all waves longer than this, 
which is typically the range covered by spectral wave models 
for sea state forecasting. 
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where 2+ = (2 — (,)/ ^J2v/a, with the surface elevation 
C,{x, t) = a cos(kx — at). 

The dissipation rate of energy is given by the mean work 
of the viscous stresses, 



p — C g p v 



(pavudu/dz) 
p w ga 2 /2 



(A6) 



Equation (A5) gives, 



du e 



yf2^ 



[cos (kx — at — z+) — sin (kx — at — 2+)] 



Figure 9. Boundary layer over waves in the absence of 
wind. Because of the larger inertia of the water compared 
to the air, most of the adjustment from the sub-surface 
velocity to the free stream velocity in the air occurs on 
the air-side of the surface. 



For a comparison with fixed bottom boundary layers, the 
leynolds number based on the orbital motion should be re- 
defined with a doubled velocity and a doubled displacement, 
i.e. Re= 4u or ba rb/^- For monochromatic waves a or b = a 
and Uorb = aa — 2na/T. For random waves, investiga- 
tions of the ocean bottom boundary layer suggest that the 
boundary layer properties are roughly equivalent to that of a 
monochromatic boundary layer defined by significant prop- 
erties (Traykovski et al. 1999). 

Although the wind was neglected here, it should influence 
the shear stresses when its vertical shear is of the order of 
the wave-induced shear. Taking a boundary layer thickness 
S and wind friction velocity u+, and assuming a logarith- 
mic wind profile, this should occur when m*/(k5) exceeds 
2u OI b/S, where k is von Karman's constant. This corre- 
sponds to, roughly, u* > u or b- For swells with T < 15 s and 



H as > 2 m (i.e. u or b,s > 0.4 m 



and winds less than 



(A7) 7ms 1 (i.e. u t < 0.2 m s x ), the wind effect on f e may 
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be small and the previous analysis is likely valid. In gen- 
eral, however, the nonlinear interaction of the wave motion 
and wind should be considered, which requires an extension 
of existing theories for the distortion of the airflow to finite 
swell amplitudes. 

Finally, the above result for the air viscosity is easily gen- 
erated to any water depth D by dividing the free stream 
velocity u+(x,t) by tanh(fcD), so that the dissipation rate 
is a factor 1/ tanh 2 (fcD) larger in intermediate water depth. 

Appendix B: Quantitative validation of H ss 

A classical analysis of SAR estimation errors is provided 
by a direct comparison of swell parameters, estimated from 
level 2 products, with buoy measurements at nearly the same 
place and time [Holt et al., 1998; Johnsen and Collard, 2004]. 

Previous validations were presented for the total wave 
height H s [Collard et al., 2005] or a truncated wave height 
H B \2 defined by chopping the spectrum at a fixed frequency 
cut-off of 1/12 Hz. For that parameter, Johnsen and Collard 
[2004] found a root mean square (RMS) difference of 0.5 m, 
when comparing SAR against buoy data, including a bias 
of 0.2 m. In the present study, we use H ss values obtained 
from both SAR and buoy spectra. 

For each wave spectrum observed in the world ocean, 
swell partitions are extracted providing estimations of H ss , 
T p , and 6 V . In practice, the L2 spectra are first smoothed 
over 3 direction bins (30° sectors) and 3 wavenumber bins, 
in order to remove multiple peaks that actually correspond 
to the same swell system. The peaks are then detected and 
the energy associated to each peak is obtained by the usual 
inverted water-catchment procedure [Gerling, 1992]. The 
swell peak period is defined as the energy-weighted aver- 
age around ±22% of the frequency with the maximum en- 
ergy. Likewise the peak direction P is defined as the energy- 
weighter direction within 30° of the peak direction. 
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Figure 10. ASAR-derived swell partition heights versus 
buoy swell partition heights after bias correction using 
eq. (3), for subset A. The solid line joins the median val- 
ues from SAR observations in each 0.1 m class of buoy- 
measured height. 



Table 1. Error statistics for swell partitions heights and 
peak periods derived from SAR wave mode data (after bias 
correction) against buoy-derived data, for subsets A, B and 
C of the co-located database. Subset A contains 2399 ob- 
servations. Subset B is restricted to UioSAR < 8m/s and 
contains 1936 observations. Subset C is further restricted to 
SAR-buoy distances less than 100 km, and contains only 460 
observations. RMSE stands for root mean square difference, 
while the NRMSE is the RMSE normalized by the root mean 
square observed value. The scatter index (S.I.) is equivalent 
to the NRMSE computed after bias removal. Finally, r is 
Pearson's linear correlation coefficient. 
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A preliminary validation of H 3S was performed by Collard 
et al. [2006], using L2 processing applied to 4 by 4 km tiles 
from narrow swath images exactly located at buoy positions. 
That study found a 0.37 m r.m.s. error. This smaller error 
was obtained in spite of a 4 times smaller image area that 
should on the contrary produce larger errors due to statisti- 
cal uncertainties. This suggests that a significant part of the 
"errors" in SAR validation studies are due to the distance 
between SAR and buoy observations. 

The swell height validation has been repeated here, us- 
ing all buoy data from 2004 to 2008, located within 200 km 
and 1 hour of the SAR observation. These co-located data 
are made publically available as part of the XCOL project 
on the CERSAT ftp server, managed by Ifremer. Because 
we wished to avoid differences due to coastal sheltering and 
shallow water effects, we restricted our choice of buoys to 
distances from the coast and the 100 m depth contour larger 
than 100 km. As a result, most selected buoys are not di- 
rectional, and we use partitions in frequency only. For the 
present validation, differently from other section of this pa- 
per, we thus define a swell partition as the region between 
two minima of the frequency spectrum. The correspond- 
ing energy E s gives the swell height H ss = A\fWl. The buoy 
swell height is then defined from the energy contained within 
the frequency band of the SAR partition. The peak period 
is then estimated as the period where the buoy spectrum 
is maximum. The database includes 15628 swell partitions 
observed by the SAR, with matched buoy swell partitions. 

Many of these observations correspond to relatively short 
swells, for which the waves are poorly imaged. We have thus 
defined a subset of the database by imposing the following 
conditions. First the image normalized variance, linked to 
the contrast intensity and homogeneity, should be in the 
range 1.05 to 1.5, which limits the dataset to 6651 obser- 
vations. This removes SAR data with non-wave features 
(slicks, ships ...) that would otherwise contaminate the wave 
spectra. Second, both the SAR and buoy peak periods are 
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restricted to the 12 to 18 s range, which reduces the dataset 
to 4136 observations, and removes most of the problems re- 
lated to the azimuthal cut-off. Third and last, the SAR- 
derived wind speed Uiosar is limited to range from 3 to 9 
m s _1 in order to remove low winds with poorly contrasted 
SAR images and high winds which may still cause some im- 
portant azimuthal cut-off and contamination of swell spec- 
tra by wind sea spectra. This gives subset A, with 2399 
observations. The resulting heights are compared to buoy 
measurements in figure 10. 

When the maximum wind is reduced to 8 m s~\ giving 
subset B, the differences between SAR and buoy data is re- 
duced, with further reductions when the maximum distance 
between SAR and buoy data is reduced from 200 to 100 km 
to give subset C (table 1). 

Appendix C: Derivation and verification 
of the asymptotic swell energy without 
dissipation 

CI. Derivation 

The swell energy E s is an integral of the local spectrum 
G over both frequencies / and arrival directions 9, 



E s (a) = 



G(t,<f>',\'j,e')due'. (ci) 



Using eq. (1) this local integral, can also be written as an 
integral over the entire source area fl. The spherical dis- 
tance between any point P(<£, A) in the source region and 
the observation point 0(a, 0) is a'. The observed frequency 
that is due to this source point is 



f = gt/(47rRa') = f ^. 

a' 

We may replace / by a' in eq. (CI), 



E s (a) = 



./'<) 



a 2 G(to,4>,A,/,fl) dg/da/ 



(C2) 



(C3) 



For a circular uniform storm of radius r with isotropic spec- 
tra, as used in figure 11, E s (a) is given by the integral over 
a' weighted by the directional width of the spectrum A9' . 
That width is given by the spherical trigonometry relation- 
ship 



A<9' = 2(9' 



2 arccos 



COS V5 max — COS a COS a 



sin a sin a' 



,(C4) 



with v5 max =r/R. 

For general spectral distribution, we may transform the 
integration variables (a! ,6') which are the colatitude and 
longitude coordinates on the sphere with a pole at the ob- 
servation point, to coordinates (ip, A) with a pole in the cen- 
ter of the swell field at to- The transformation Jacobian 
is simply given by the equality of the elementary area on 
the unit radius sphere dA — sin^d^dA = |cos0|d0dA = 
sina'da'd#'. We thus have 



E.(a) 



_J^_ / a sin a Q ( Q) ^ ^ 

asm a _/ n a * sin a 



(C5) 



G(to,<j>,\,fo,9o)dA 



l + O 



+ f \G(to,4>,\,f,9)-G(to,4>,\,fo,6)]&A 



+ J [G (to, </>, A, f ,0) - G (to, <p, A, f , 0o)] dA j . 

(C6) 

where 9 is the direction of the great circle at the genera- 
tion point that goes through that point and the observation 
point. 9 is thus a function of <f>, A, a and 9o- A a is maxi- 
mum value of \a' — a\, i.e. the radius of the source region 
divided by the Earth radius. 

For continuous spectra, the second integral on the right 
hand side of eq. (CI) goes to zero as a goes to infinity 
(which, on the sphere is limited by it) since the part of the 
source spectra that contribute to the observations shrink to 
a smaller and smaller neighborhood around fo and 9q. The 
observed frequencies / are limited by 



l/-/o| < Ao 
fo a 



(C7) 



This is enough to guarantee that this second integral also 
contributes a deviation £2 from the asymptote, limited by 



. , A a dG 

e 2 < A / max — 

a I df 



(C8) 



where the maximum is taken over all the contributing com- 
ponents. 

Similarly, the outgoing directions 9 received at the obser- 
vation point are also limited to a narrow window as a in- 
creases, giving another deviation term £3. Using the sine for- 
mula in the triangle OPN, sin(?/sin(7r — a) = sin A/ sin a'. 
Thus, in the far field of the storm and its antipode, 9 is close 
to 6*o = A. Thus sin 9 — sin do = sin A (sin a — sin a') / sin a' , 
which is less than A a / sin a', and therefore \0 — 9o\ is less 
than arcsin (A a / sin a'). If one does not get too close to the 
storm or its antipode (say, A a < a < tt — 2A a ) then we 
can give an upper bound of 1/sino/ as a function of a and 
obtain 



£3 < A arcsin (2 ^ a \ 
\ sin a ) 



dG) 

!,iax l¥/' 



(C9) 



The deviation from the asymptote due to £3 is thus of the 
order of A^/sina and may increase close to the antipode, if 
waves from a wide range of directions can reach that point. 
In practice this does not happen since continents and is- 
land chains block most of the arrival directions at the an- 
tipode, leaving only a small window of possible arrival di- 
rections [e.g. Munk et al, 1963]. Directional wave spectra in 
active generation areas are generally relatively broad with 
dG/d6/G typically less than 2 for directions within 30° of 
the main wave direction. On the contrary, typical storm 
spectra can give fdG/df/G as large as 10 for frequencies 
within 30% of the peak frequency. We thus expect the de- 
viation £2 to dominate over £3. 

C2. Verification 

For a useful comparison with observations, the asymp- 
totic swell height evolution should be approached on a scale 
smaller than the ocean basin scale. This is easily tested 
for storms with spatially uniform spectra over a radius r, 
by evaluating the integral in eq. (C3). We chose a cen- 
ter frequency fo and consider the swell energy at a distance 
Ra + A x and a time t(a) such that Ra = gt/(4ivfo). A x 
is thus an error relative to the theoretical position of the 
wave group of frequency fo- The relative difference of E s (a) 
and its asymptotic evolution l/(a sin a) depends only on the 
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Figure 11. Convergence of the swell energy E s integral 
(CI) towards the asymptote l/(asina), as a function of 
spectrum width for a storm diameter of 1000 km. The 
result is independent of the choice of f p . In practice the 
calculations were made for f p = 0.07 Hz (T p = 13 s). The 
three lines for each case correspond to position errors A x 
of -200, and 200 km relative to the great circle trajec- 
tory. For all cases considered here the deviation is less 
than 20% beyond 4000 km from the storm center. 



spectral shape in the storm, the relative frequency fo/f P 
where f p is the peak frequency, the distance a, the storm 
size a/(r/R), and the position error A x . Results are shown 
in figure 11 for isotropic directional spectra, in which case 
£3 = 0. Although an isotropic spectrum is not realistic at 
all, it allows for simple calculations. As discussed below, the 
effect of directional spreading is expected to be less impor- 
tant than the shape of the frequency spectrum. 

As indicated by eq. (C8), the contribution of both the 
spectral shape and fo / f p comes through the maximum rela- 
tive variation of F in the frequency interval that contribute 
to E s . Here we take the spectra in the storm to have a 
JONSWAP shape [Hasselmann et al, 1973], that we adjust 
by varying the peak enhancement factor 7. A relatively 
broad Pierson and Moskowitz [1964] spectrum is obtained 
with 7 = 1. If we chose fo = f p , this spectrum will give 
smaller deviations of E s from the asymptote (solid lines in 
figure 11), than narrower spectra with larger values of 7. 
Young [2005] showed that wave spectra in Hurricanes gen- 
erally fall in between these two categories, with a typical 
value 7,7 = 1.7. Larger deviations from the asymptote are 
obtained for fo < f p , since the forward face of the spectrum 
is very steep, while smaller errors are obtained for fo > f p 
due to the more gentle decrease of F towards the high fre- 
quencies. 

Similarly, large deviations are produced if the observa- 
tions are made in a direction far from the peak generation 
direction in cases when the directional spectrum is narrow. 
For observation directions 30° from the peak direction, and 
spectra with a cos 4 directional distribution, the deviations 
are still dominated by the dispersive term £2, as expected. 

The other important factor is the distance a relative to 
the storm size (pmax = r/R. A faster convergence is ob- 
tained for smaller storms. If observations do not correspond 
exactly to a theoretical propagation at a group speed C g but 
are within a distance of the theoretical position, the val- 
ues of -Es will also be affected in a way that depends on the 
spectral shape. In the absence of energy gains or losses, and 
for realistic storm sizes and spectral shapes, the deviation 
of observations from the asymptote should be less than 20% 
for x > 4000 km when A x < 200 km is enforced. 
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